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Neutrino cloud instabilities just above the neutrino sphere of a supernova 
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Most treatments of neutrino flavor-evolution, above a surface of the last scattering, take identical 
angular distributions on this surface for the different initial (unmixed) flavors, and for particles and 
antiparticles. Differences in these distributions must be present, as a result of the species-dependent 
scattering cross-sections lower in the star. These lead to a new set of non-linear equations, unstable 
even at the initial surface with respect to perturbations that break all-over spherical-symmetry. 
There could be important consequences for explosion dynamics as well as for the neutrino pulse in 
the outer regions. 


The study of collective neutrino interactions in the out¬ 
skirts of the supernova explosion has generated an explo¬ 
sion of publications These collective effects come 

about as a consequence of the standard model neutrino- 
neutrino coupling, and change the correlations between 
neutrino flavor and neutrino energy. In doing so they 
may profoundly affect the neutron-proton ratio, and R- 
process rates in layers far above that of the last scatter¬ 
ing. 

In these studies what is added to the observationally 
well established single v oscillation and matter effects can 
be described as totally forward, standard model, interac¬ 
tions among the v, viz., v{q)-\-v{k) iy{q) + iz(k). Noth¬ 
ing except iz-flavor enters into the process. In principle, 
therefore, one can follow outward each neutrino, contin¬ 
uously updating the flavor density-matrix of each to take 
into account the interactions with all of the others. Con¬ 
sideration of this process, as codified in a Liouville-type 
equation [ij , has led to a number of large-scale numerical 
calculations, e.g. refs. i-0, in which the underlying 
assumptions have been named the ‘T-bulb” model. 

These calculations begin at the surface of a neutrino- 
sphere, a loosely-defined surface beyond which scattering 
is considered ignorable, and they assume spherically sym¬ 
metrical neutrino flows from below. On this i^-surface 
the v angular distributions are taken as uniform, with no 
variation in angle over the outward directions, but with 
no inward bound v^s. The initial states for these neutri¬ 
nos are taken as “flavor-diagonal”, that is with no mixing 
of Ve, or of Here is that combination and 

Vt that is best suited to two-neutrino simulations. 

There have been a number of recent works [13 -[33 de¬ 
voted to the possibility of “azimuthal” or “multi-angle” 
instabilities in the above picture; that is, that some small 
perturbation of the v flavor density matrix that breaks 
the spherical symmetry of system can grow exponentially 
beyond some radius of onset. In these examples the in¬ 
stabilities would occur at distances of 100 km. or more 
from the ‘V-surface.”, a sphere of last scattering that is 
at a radius of about 15 km.. By contrast we find in¬ 
stability at or very near this surface. There is a reason 
for this difference: In some imprecise way it is already 
known that the onset of a “multi-angle” instability, stem¬ 


ming from the v — v interaction, requires a state in which 
angular distributions are flavor dependent [ 13 , 0 , and 
in the references cited above these develop through joint 
action of the v — v interaction and neutrino oscillation 
terms. The latter are greatly suppressed by the interac¬ 
tion with electrons until the z^’s are out to the 100 km. 
region. In contrast we begin at the z/-surface with a phys¬ 
ically based flavor dependence in the angular distribution 
that makes the system unstable to the tiniest seed from 
oscillation effects. 

Some features of our approach are: 

1. We follow custom by taking distributions on initial 
surfaces to be angle-independent in the outgoing hemi¬ 
sphere. But we take the radii of these surfaces to depend 
on “flavor”, which in our terminology will distinguish 
{i'e,9e,Vx,iix}. This is motivated by the fact that the 
average last-scattering radius r ,^^, for z/g, is considerably 
larger than that for Pg which, in turn, should be larger 
than that for zzj,, ■ 

2. Then we begin our real calculation (still with flavor- 

diagonal distributions) at an initial distance The 

other three favors on this surface have narrower angular 
cones than does z/g, just from the effects of free streaming 
from their lower surfaces of origin. 

3. The distributions on an initial surface that are pro¬ 

duced from the above are perfectly spherically symmetric 
in the sense that zz angular distributions, f{0), where 9 is 
the angle with respect to the local normal, are the same at 
all points of the sphere. Testing the non-linear evolution 
equations for stability against non-spherically symmetric 
perturbations in the flavor density matrices, we find in¬ 
stability at the initial surface at tq (with rp > ) for an 

interesting range of parameter values. 

4. The system described above still has perfect spher¬ 
ical symmetry in its flavor distribution and flow. The 
perturbations that grow, at a rate of order Gpni/ as it 
turns out, are characterized by having a combination of 
tiny flavor mixing and spherical asymmetry. Seeds of 
such irregularities can come about in many ways. We 
mention one that does not even require asymmetry in 
the v number flow: Neutrino oscillation terms will give 
tiny flavor mixings during the propagation between last 
scattering and the starting surface of the calculation. Ir- 
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regularities in electron densities in this region will lead 
to angular irregularities in the starting flavor-density ma¬ 
trix, insignificant in stable cases but crucial in unstable 
cases. 

We begin from the Liouville equation 3l| , 32| for time 
independent flow, 


n • V(7Q(r, n) = 2^/^(27r) x 

f dfln>a, 3 (r, n')(l - n • ri') + {O&M} , (1) 

where the vectors indicated by arrows and the vector 
product are in i/-flavor space and fi is neutrino velocity 
in real space. {O&M} stands for the standard oscilla¬ 
tion and matter terms, linear functions of the compo¬ 
nents of da , which we drop from the development below, 
since they will enter only in the preparation of the ini¬ 
tial state for the calculation. The index a keeps track 
of distinguishing the origin of a neutrino as a pure fla¬ 
vor state in one of the groups {i/g, i/g, Pa,}. The set 
Ca = puts in minus signs for anti-z/’s inside 

the sum. In ([T|), and everything that follows, the energy 
spectrum is irrelevant; it is number densities that enter. 

In the case of a spherically symmetric system we have 
(Ta(r, fi) = dair,r ■ fi), whence, 


o'o,a = fa{r sin 0), where the /q’s are any functions. The 
“neutrino bulb” model assumes (To,a = Ua0(ro — rsind), 
independent of a, where 0 is the Heaviside function and 
na are the number densities of the species a on a common 
surface. We now relax this assumption by taking differ¬ 
ent effective neutrino-surfaces for the different groups of 
neutrinos. 


o'o.a = na&{ra — T sin 6 ), 


( 5 ) 


where }i with being the largest. 

Next we re-express the right hand side of m in a series 
of steps: 

a) . We replace cos 9 by a new variable s with the re¬ 
lation chosen differently for each of the four functions 
<i>Q,(r, COS0), where cos0 = 1 — sha{r), and 

/i„(r) = l-[l-(r„/r)2]i/2. (6) 

b) Step a) above entails the replacement in (g]), 

J c?Hft/cro,/ 3 (r, cos 0 ')$/ 3 (r, cos0')(l — n ■ h') ^ 

2-Knp Jq ds'$p{r, s')hp{haS + hps' — hahpss') (7) 

c) . In the new variables we have. 


n- 


[n-V (?-,?• n) = 

dida (r, f • ' n]) 92 da (r, f • fi^ 

da (r, cos 9) = Oda {r, cos 9) 


cos 9— + ^ ^_ 

cos (7 -I- ,, a cos e 


( 2 ) 


where di and 92 are the derivatives with respect to the 
first and second arguments of da, the angle 9 is between 
r and n, and the last line defines the operator O. Using 
the r dependent variable of would have eliminated 
the 9/9 cos 0 term in O, but made the p — p interaction 
very difficult to deal with when the P'-sphere radii depend 
on flavor. Results the equivalent of ([2]) are found in the 
appendix to ref.ls^. The da’s are related to the flavor- 
density matrices, cos 9), for individual neutrinos by 


da{r, cos 9) — ao^a{r,cos9) ^a{r,cos9), (3) 

where ao^a{r,cos9) is the total neutrino density in the 
group a. We have Oda{r, cos 9) = 0 in the absence of iz 
interactions, and OuQ^air, cos9) = 0 even in the presence 
of the neutrino terms, so that o becomes, 

0$Q(r,cosd) = 2^/^(27r)“^Gi7’$a(7',COS0) x 

^ E/ 3 C /3 / dHn/CTo,/3(7’, cos6»')^/3(i’, n')(l - n • n'). (4) 

In the flavor diagonal initial state we have = 1, 
^ (Here and in what follows the 
superscripts x,y,z stand for directions 1,2,3 in flavor 
space.) The solutions to Oao^a(j,n) = 0 are of the form 


0 = [l-sha{r)]£ + 


d r 

{s'^hoL {r) — s)h'^{r) 

M7) 


ha{r)s^ —2s 


ds ’ 


and the evolution equation is, 

C><i)„(r, s) = 

s,) x Ua{r, s ), 

where 

Ua{r,s) = fg ds'lska + s'hp - 

ss' hahp\^p{r, s'). 


( 8 ) 


(9) 


( 10 ) 


We shall examine instabilities by making the perturba¬ 
tion $ — $ -I- X , then linearizing the RHS of ([9]) in y. 
In all calculations we shall discretize in the variables s, s' 
with N equal steps in the interval (0,1). 

This calculation will be local, in the sense that the 
instability criteria will involve the background flavor- 
density matrix only at a single point, and for that we 
use the results of free-streaming from the several sepa¬ 
rate z/-surfaces. When we then look to plot the growth of 
instabilities, we concern ourselves with a volume of sur¬ 
rounding space small enough such that in the absence of 
instabilities there would have been no appreciable varia¬ 
tions of $. From ([9]) we obtain. 


OXa,i 


j- 3/2 


Gf 




XoL,i X X],a=l = l ^ 0,j^OL,i-,PJ 


X J2g=l Sfcl X/3.: 


-b 

( 11 ) 
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where Wa,i;/ 3 j, with discretized angular indices (i,j), is where 
given by 




hi3ni3C/3[ihaN '^ + j hpN ^+ijhahpN 

( 12 ) 


and 


,(i) 




P',3' 


(17) 


In all that follows, we take = 1, = 0 in (fTTl) 

for all {a,*}, since the initial state is to have seen no 
neutrino oscillation, giving. 


where. 


Ox‘S = -Y,«:,-m4: 


(x) 

3 


= 2 - 3 / 27 r- 2 Gi. 7 V-i x 

^ct,P ^i,j ^o,,'^‘,P', 3 ' * 1 / 3,1 


(13) 


(14) 


The indices{a,i} and {P,j} henceforth will be consol¬ 
idated into a and &, respectively, in a 4iV dimensional 
space. The eigenvalues of Ma^, as a matrix in this space, 
then indicate instability if any have an imaginary part. 
(It is the X —y symmetry in (1131) that converts an 
imaginary part here effectively into a real part for an 
eigenvalue of the system as a whole.) The coupled equa¬ 
tions dni) would then have exponentially increasing solu¬ 
tions, as long as an initial value has any component in the 
direction of the eigenvector. The first term on the RHS 
of (03 is diagonal and real but the second term is non- 
hermitian; nonetheless for a wide range of parameters 
and ria that we have considered we find stability. We be¬ 
lieve this stability with respect to spherically symmetric 
perturbations to be a general result. 

To address azimuthal instabilities at the initial sur¬ 
face r = we take the unperturbed ^a,i as given 
by our previous expressions, and therefore with angu¬ 
lar distributions that are functions of as above. 

For the linear azimuthal perturbation we take the form 
= ;yf“’i'l[cos6()] cos(/>, where the azimuthal an¬ 
gle (j) is defined locally with respect to a tangent to the 
surface. Now in the analogue to the final term in (HU, 
performing the 4>’ integral over fi ■ n', we obtain, 


Ox':l = [xi»> ■I-S 

E/j j X% & "11 '■S "tj 




(15) 


where the discretized sin 9 factors from n • n' have the 
form, Va,i{r) = [1 — (1 — iha{r)/N)'^]^/'^. Following the 
previous development we then recapture (USD, but with 
M replaced by. 




.3=2 


-3/2,,-2 


GfN 


-ir 


l,j + '"a.i;/: 


(16) 


=2 Hpnph^ Vc,,i vpj . (18) 

We shall look for instabilities at r = , the outermost 

neutrino-surface. Their presence now depends on the val¬ 
ues of four parameter ratios, riv^lnu^, 

JT-j/j,/fz-i/e ■ Explosion calculations give rather definite val¬ 
ues for the number density ratios, which can be calcu¬ 
lated from the ratios of luminosity to average energy. The 
raw ratios use the luminosities given by the simulations 
at a common radius outside all neutrino spheres, but still 
in the region in which neutrino oscillations matter very 
little. But the are the densities at the respective v- 
surfaces; and we must therefore apply a correction factor 
to the raw ratios for a = Vf.,Vx,i>x- 

Here we shall use the values given in fig. 7 of ref. in 
the range of times, .3 sec.-l sec., post-bounce. We find 

nujni,^ Ki .77[r^,/rpj2^ Ki .62[r^^/r^J2_ 

have not found as authoritative predictions for neutrino- 
surface radii; a neutrino-surface is an idealization that is 
not contained in the data. But we believe that, given the 
idealization, the opacities and the matter density profile 
dictate that Within this range and 

for the above values of ria ratios, we find instability in a 
region approximately described as, 

> .U+ .bb Jr. (19) 

Thus if the radius for is 80% of the radius for z/g then 
only the region > .88 is unstable. 

The instability is manifested by imaginary parts in a 
single physical pair of eigenvalues. A, A*, of M converging 
rapidly to a limit as we increase the subdivision number 
N (say, from 16 to 128), while spurious imaginary parts in 
other eigenvalues as a result of finite subdivision (noted in 
[ 3 ^) go to zero as As an example, when we choose 
?T-Pe = 10^^cm“3 at r = 15 km, and choose the ratios 
rp^ = .93 Tp^, rp^ = .8rp^ we find Im[A] = ±.008G_Fnp^ ~ 
.32 (meter)“^. 

The critical region given by m appears to be likely 
to be realized at some times or places during the event. 
As to “places”, we suggest a broader context for our con¬ 
siderations. The interior simulations suggest an environ¬ 
ment that is locally very irregular in its parameters of 
density and temperature. This will translate into very 
irregular neutrino angular distributions as well. If, as 
in simulations like those of i-0, the important flavor 
evolution takes place over a region in r many times the 
size of the z/-surface, these irregularities may get aver¬ 
aged out, since at most points on the way out we see 
solid angle of nearly 27r of the zz-surface. But if the scale 
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of V flavor exchange is tens of meters instead, then we 
see only a small fraction of the complete (/-surface as we 
sit only a little ways above it, and the averaging argu¬ 
ment is less sustainable. We regard our present work as 
a step of interpolation between the methods and results 
for the usual smoothed-out model and those of a future 
statistical model. 

For the remainder of the paper we shift to taking the 
system at time t = 0 to have been established in a state 
of steady symmetric flow, and we then set O — dfdt, 
turn on the perturbation and follow the time develop¬ 
ment at a fixed point, seeking the results for some pe¬ 
riod of time short compared to rg but long compared 
to the exponential rise times. Following that path we 
now examine the back-reaction on of the growing 
transverse modes, . yve use the complex repre- 

sentation, Xa = Xa + iXa in what follows. Taking the 
normalized eigenvector of M for the growing mode to be 
^a, with eigenvalue A, we assume Xa{t) ~ for all 

a = that is, we take the growing mode to domi¬ 

nate the (a;, y) components. The evolution equation (U) 
couples all of the ' modes to the perturbations, y. 
Defining z(t) = 4)1^^ (i. e. summed over all 

angles and all flavors), we find, 

z{t) = 

= |g(t)pIm[A], (20) 

where to get the second form we hrst replaced by 
since is Hermitian, then used the fact 
that is an eigenvector of -I- with A the grow¬ 
ing mode eigenvalue of the matrix M of (fTHl) . This is 
supplemented with, 

g{t) = iXz{t)g{t) (21) 

which comes directly from m- Using (1^ we obtain the 
first equality in 

|[| 5 (t)n = -2Im[A]z(t)|ff(t)p = -2z{t)m (22) 

the second following from (1^ . and giving the solution 
Iff (OP = + C”. We take z(0) = 1, so that C = 

|g(0)p -I- 1. Now differentiating (20) and using the above 
we obtain, 

m = 2[iuv\]Mtf - (1 + \mfMt) . ( 23 ) 

The parameter |g(0)p is the square of the initial v mix¬ 
ing amplitudes, from whatever source below in the star, 
in particular the one discussed in 4) in the introduc¬ 
tion. In fig. 1 we plot the results for values |(7(0)p = 
10“^, 10“®, 10“®, using the eigenvalue calculated previ¬ 
ously to express the scale in meters/c. 

Recall that in our definition, ) tracks the evolution 
in the stream that has original flavor a. A transition in 


z 

1.0 

0.5 


-0.5 
- 1.0 

FIG. 1: 

Flavor turnovers as measured by z(t), plotted against time. 
The dotted, solid, and dashed curves are for the cases 
lff(0)P = 10“"^, 10“®, 10“®, respectively. 



which (4Af)“^ 'Yua goes from 1 to —I as in hg.l means 
complete e •<->■ x exchange for each beam. The curve dur¬ 
ing the short transition is almost a perfect “tanh” kink 
as long as the initial values for |g(0)p are sufficiently 
small that the length of the initial plateau is long com¬ 
pared with the transition time, as in the plotted cases. 
Note that the plots begin at the middle of a plateau and 
the turn-over time is proportional to |(log |ff(0)|)|“^. For 
each reduction of a factor of 10 in the size of the seed, 
( 7 ( 0 ), the turn-over is 2.5 meters farther away from the 
initial surface. This observation is important in dealing 
with the question of the matter effects, which in other 
calculations damp the growth of oscillations in the re¬ 
gions of higher electron density. In our calculation they 
barely matter; a tiny push, |( 7 ( 0 )p, gets the unstable sys¬ 
tem moving, but damping this push by a factor of 100 
hardly changes the results. The non-linear interaction 
exchanges the directions and energies of e and x neutri¬ 
nos without changing their numbers, and likewise for the 
i/’s. The total interaction energy with electrons in the 
medium therefore does not change in this transformation 
since it depends just on densities and not on spectra. 

This is quite different from the situations in other liter¬ 
ature in which one starts on a single surface in a flavor- 
diagonal state with a species independent angular dis¬ 
tribution and finds instabilities only after a substantial 
amount of mixing from the oscillation terms. At the 
shorter distances (say 15 km) with high electron densi¬ 
ties this mixing is strongly suppressed by the interaction 
with electrons. This is why in typical calculations one 
gets sudden turnovers only at distances of order 100 km. 

Fig. 1 shows total (/g O interchange over a very 
short time from the beginning of the calculation (for a 
case where Ve = Vx). Since the temperature ratios at the 
(/-surface are expected to be of the order ~ 2 

this would represent a 100% heating for the j/g, which in 
turn would much more efficient at heating the matter in 
this region. This sudden heating in the neutrino-surface 
region would be expected to be very important to the 
explosion calculation. We can only speculate about what 
might happen in a more realistic case in which the gross 
governing feature remains the differing neutrino-surface 
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radii among the species, but the seeds (still on a nearly 
planar local initial surface) are more irregular on that 
surface, e. g. in ^ dependence. As we move outward the 
different domains originating with different areas on the 
original surface will begin to merge. All of this appears 
most likely to lead to a rather quick averaging in which 
we settle at z ~ 0 on average in fig. 1. Then instead of 
a 100 % gain in temperature of t'e ’s we would have in 
effect a 50% increase, still very potent in terms of heating 
the surrounding matter. 

We believe that this last possibility is sufficiently likely 
to justify a complete explosion calculation in which the 
flavor dependence of neutrino distributions is homoge¬ 
nized computationally, with respect to flavor-energy dis¬ 
tributions, at frequent intervals, as the calculation pro¬ 
ceeds outward through the neutrino-sphere region. The 
object would be to test qualitatively whether or not an 
extension of the considerations of the present paper can 
have an important effect on explosion dynamics 

It is interesting that the structure of the equations that 
enter the above calculations is very similar to those for 
non-linearly interacting photons, where now flavor is to 
be replaced by photon polarization in real space. This 
is another reason to study general behaviors and compu¬ 
tational algorithmic development. Possible applications 
range from laboratory studies in media [s^ to the polar¬ 
ization signal in the prompt emission from gamma ray 
bursts. In the latter case the vacuum 7 — 7 interaction 
remains effectively strong out to regions of small optical 
depth, so that the data could be useful in understanding 
the source. 

I am much indebted to Georg Raffelt for making a key 
observation at the very beginning of this work. 
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